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Abstract 

In  this  research,  a  Simulink  model  of  a  standalone  vehicular  solid-oxide  fuel  cell  (SOFC)  auxiliary  power  unit  (APU)  is  developed.  The  SOFC 
APU  model  consists  of  three  major  components:  a  controller  model;  a  power  electronics  system  model;  and  an  SOFC  plant  model,  including 
an  SOFC  stack  module,  two  heat  exchanger  modules,  and  a  combustor  module.  This  paper  discusses  the  development  of  the  nonlinear  dynamic 
models  for  the  SOFC  stacks,  the  heat  exchangers  and  the  combustors.  When  coupling  with  a  controller  model  and  a  power  electronic  circuit  model, 
the  developed  SOFC  plant  model  is  able  to  model  the  thermal  dynamics  and  the  electrochemical  dynamics  inside  the  SOFC  APU  components,  as 
well  as  the  transient  responses  to  the  electric  loading  changes.  It  has  been  shown  that  having  such  a  model  for  the  SOFC  APU  will  help  design 
engineers  to  adjust  design  parameters  to  optimize  the  performance.  The  modeling  results  of  the  SOFC  APU  heat-up  stage  and  the  output  voltage 
response  to  a  sudden  load  change  are  presented  in  this  paper.  The  fuel  flow  regulation  based  on  fuel  utilization  is  also  briefly  discussed. 

©  2006  Elsevier  B.Y.  All  rights  reserved. 
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1.  Introduction 

The  solid-oxide  fuel  cell  (SOFC)  is  a  type  of  fuel  cell  in  which 
the  electrolyte  is  a  solid,  nonporous  metal  oxide.  The  SOFC 
operates  at  very  high  temperatures  (600-1000  °C).  In  the  cells, 
oxygen  ions  are  transferred  through  a  solid-oxide  electrolyte 
material  to  react  with  hydrogen  on  the  anode  side.  Their  off¬ 
gases  are  often  used  to  fire  a  secondary  gas  turbine  to  improve 
electrical  efficiency. 

The  SOFC  auxiliary  power  unit  (APU)  draws  the  interest 
of  many  researchers  because  of  its  high  efficiency,  low  emis¬ 
sion,  modular  structure,  and  high  generation  capability.  Many 
SOFC  applications  are  stationary.  Recently,  the  development  of 
an  SOFC  APU  for  vehicular  usage  on  long-haul  trucks  is  gaining 
more  and  more  attention.  In  this  paper,  we  modeled  and  simu¬ 
lated  a  completely  standalone  SOFC  unit  fueled  by  hydrocarbon 
fuel  for  control  applications. 

The  technical  challenges  of  designing  an  SOFC  APU  con¬ 
troller  are 

•  Thermal  management. 
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•  Fuel  management. 

•  Transient  and  load  following  characteristics. 

Thermal  management  is  critical  in  maintaining  the  normal  oper¬ 
ation  of  an  SOFC  system.  SOFCs  operate  at  temperatures  rang¬ 
ing  from  600  tol000°C.  A  standalone  SOFC  unit  should  be 
able  to  run  with  fuel  and  air  supplied  at  their  normal  storage 
temperature.  The  high- temperature  working  condition  requires 
that  the  inflow  gases  need  to  be  heated  to  their  working  tem¬ 
perature  prior  to  being  fed  to  the  SOFC  stacks.  As  a  standalone 
unit,  the  SOFC  system  should  be  able  to  heat  the  inflow  gases 
using  the  energy  generated  by  combustion  of  a  portion  of  the 
inflow  gases  and/or  the  exhaust  gases.  Appropriate  amounts  of 
fuel  and  air  heated  to  their  working  temperature  should  be  guar¬ 
anteed  without  any  extra  facilities,  other  than  the  SOFC  unit 
itself. 

Fuel  management  allows  the  SOFC  to  run  at  a  higher  effi¬ 
ciency,  while  the  voltage  is  maintained  within  the  desired  range. 
The  transient  and  load  following  characteristics  of  an  SOFC 
stack  determine  the  design  of  the  power  electronic  circuits. 

One  of  the  models  discussed  in  this  paper  is  the  electrochem¬ 
ical  dynamic  model.  The  electrochemical  reaction  happens  in 
the  SOFC  stack.  The  electrical  output  of  the  SOFC  can  be  influ¬ 
enced  by  power  load,  fuel  and  oxidant  gas  flow  rate,  temperature, 
and  partial  pressure  of  the  reacting  gases.  The  electrochemical 
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Nomenclature 

A  area  (m2) 

Ci  specific  heat  (kJ  kg  - 1  K- 1 ) 

C  discharge  coefficient 

Ci  thermal  capacity  (J  K-  1 ) 

7s°  standard  reversible  cell  potential  (V) 

F  Faraday’s  constant  (96485.3415  C  mol-1) 

G  Gibb’s  free  energy  (J) 

AG  change  of  Gibb’s  free  energy  (J  s-1) 

77  enthalpy  (J  mol  - 1 ) 

AH  total  enthalpy  change  caused  by  chemical  reac¬ 
tions  (Js-1) 

I  current  (A) 

K  specific  heat  ratio  of  the  discharged  gas 

Kv  equilibrium  constant  of  the  reaction 
mi  mass  of  the  flow  or  one  SOFC  cell  (kg) 

M  gas  molecular  weight  (g  mol  - 1 ) 

n  number  of  moles;  or,  number  of  electrons  gener¬ 

ated  in  the  cell  reactions 
rti  molar  flow  rate  of  species  i  (mol  s-1) 
n]  total  rate  of  production  of  species  i  (mol  s-1) 
ftstack  number  of  cells  in  a  stack 
Nq  number  of  cells  in  stack 

Pi  partial  pressure  of  each  component  (species  i)  (Pa) 

P  absolute  pressure  (Pa) 

Pd  downstream  ambient  (Pa) 

Pe  power  consumed  by  the  electrical  loads  (J  s-1) 

Ps  absolute  source  pressure  (Pa) 

q  heat  flow  (Js-1) 

Q  heat  flow  rate  generated  by  chemical  reactions 
(J  s-1);  or,  electric  charge  (C) 
r  ohmic  resistance  (£2) 

R  gas  constant  (8.31451  kJkmol-1  K-1) 

Ri  thermal  resistances  (K  s  J- 1 ) 

S  entropy  (JK-1) 

T  temperature  (K) 

V  compartment  volume  (m3) 

Fanode  over-potential  caused  by  diffusion  in  the  cell 
anode  (V) 

Vb_v  Butler- Volmer  over-potential  (V) 

Kathode  over-potentials  caused  by  diffusion  in  the  cell 
cathode  (V) 

Vi  volume  of  the  flow  (m3) 

Kfernst  open-circuit  Nemst  potential  (V) 

Vr  voltage  drop  caused  by  the  material  resistances 
(V) 

Wi  mass  flow  rate  (kg  s-1) 

W  electric  work  (J  s-1) 

y  -n  inlet  flow  molar  fraction  of  species  i 

Greek  letters 

oti  stoichiometric  coefficients  of  species  i 

s  heat  transfer  rate 

Pi  density  of  the  flow  (kg  m-3) 


dynamic  model  addresses  all  these  factors,  which  influence  the 
electrical  output. 

There  is  a  great  deal  of  literature  addressing  individual 
aspects  of  this  modeling  problem,  such  as  the  thermo-fluid  mod¬ 
eling  [1],  flow  modeling  [2],  and  the  modeling  of  dynamic 
responses  of  fuel  cells  [3-9].  However,  very  few  researchers 
have  developed  a  complete  system  model  for  a  standalone  SOFC 
stack,  including  dynamic  models  for  auxiliary  equipment  such 
as  the  combustor,  heat  exchangers,  the  reformer  and  power  elec¬ 
tronic  circuits.  A  complete  system  model  allows  researchers  to 
test  different  control  schemes,  adjust  design  parameters  of  each 
component,  observe  dynamic  responses  of  the  APU,  and  tackle 
the  three  technical  challenges  in  controller  design  mentioned 
above.  Therefore,  it  is  critical  to  build  such  a  model  for  the  pur¬ 
pose  of  the  SOFC  APU  control  applications. 

This  paper  is  the  first  of  a  series  of  articles  on  this  topic, 
in  which  we  present  a  Simulink  model  comprising  four  mod¬ 
ules:  the  SOFC  stack  module,  the  combustor  module,  the  heat 
exchanger  module,  and  the  controller  logic  layer,  as  shown  in 
Fig.  1.  It  is  assumed  in  this  paper  that  the  SOFC  stack  is  fueled 
by  hydrogen,  so  that  no  reformer  is  included  in  the  model.  The 
focus  of  the  paper  is  to  derive  mathematical  models  for  each  indi¬ 
vidual  component  and  integrate  these  individual  models  into  a 
complete  SOFC  model  using  Simulink.  A  control  logic  system 
as  a  separate  layer  in  the  Simulink  model  was  developed  to  mon¬ 
itor  the  normal  working  condition  of  the  SOFC  stacks,  regulate 
the  gas  flows  according  to  power  consumption  requirement,  and 
control  the  combustor  and  the  heat  exchangers.  The  simulation 
of  the  power  electronic  system  load  following  characteristics 
will  be  presented  in  follow-up  publications,  in  which  different 
controller  design  schemes  and  the  module  of  power  electronic 
circuits  will  be  presented. 

2.  Calculations  and  results 

In  Sections  2. 1-2.3,  we  present  the  SOFC  stack,  the  com¬ 
bustor  and  the  heat  exchanger  modules.  Simulation  results  are 
reported  in  Section  2.4. 

2.7.  The  SOFC  stack  module 

Two  dynamic  responses  are  considered  in  modeling  the 
SOFC  stack:  the  thermal  dynamics  and  the  electrochemical 
dynamics. 

2.7.7.  Modeling  thermal  dynamics  inside  an  SOFC 
stack 

To  calculate  the  SOFC  stack  temperature,  a  detailed  and  accu¬ 
rate  computational  fluid  dynamics  (CFD)  model,  as  shown  in 
Fig.  2(a),  can  be  used  to  simulate  every  cell  in  great  detail  and 
to  obtain  the  temperature  gradient  of  each  cell  [1].  However, 
this  method  is  quite  involved  and  takes  a  significant  amount  of 
computational  time  to  complete;  thus,  it  is  not  applicable  for 
dynamic  modeling  of  an  SOFC  system  and  real-time  controller 
design.  For  control  applications,  it  is  not  necessary  to  have  such 
great  accuracy  because  the  feedback  will  correct  a  considerable 
amount  of  error  in  the  model.  Therefore,  it  makes  sense  to  com- 
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Fig.  1.  The  block  diagram  of  a  solid-oxide  fuel  cell  auxiliary  power  unit. 
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Fig.  2.  (a)  Fuel  cell  stacks,  (b)  a  distributed  thermal  model,  and  (c)  a  lumped  thermal  model. 
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bine  the  spatially  varying  temperatures  into  single  temperatures 
Tair,  7fuei  and  rsofc  for  the  inflow  air,  fuel,  and  the  SOFC  stack 
temperatures,  respectively.  Another  temperature  that  was  used 
in  the  SOFC  modeling  is  7internal>  which  represents  the  temper¬ 
ature  inside  a  cell,  where  chemical  reactions  happen.  It  is  noted 
that  temperature  sensors  cannot  measure  7internai  directly.  Csofc 
is  the  sensor-measured  SOFC  stack  temperature  signal,  which 
can  be  sent  to  the  controller.  The  simplification  results  in  a  dis¬ 
tributed  lumped  model,  as  shown  in  Fig.  2(b).  An  aggregated 
lumped  model  can  be  obtained  by  aggregating  the  distributed 
lumped  model  for  each  cell,  as  shown  in  Fig.  2(c). 

The  heat  flow  q  transferred  from  the  air  and  fuel  to  the  SOFC 
stack  can  be  represented  as 


<7  air  —  WJair^va(Fair 
^fuel  =  ^fuel  ^'vf  ( ^fuel 


1 

^Internal)  =  (^air  —  internal)? 

R\ 

1 

—  ^Internal)  —  (Cfuei  —  ^Internal) 

Rl 


where  R\  = 


l 


and  R2  = 


1 


^air^va  ^fuel^vf 

Also,  u;air  and  u;fuei  are  the  air  and  fuel  flow  rates,  respec¬ 
tively,  cva  and  cvf  the  specific  heat  of  air  and  the  fuel,  respectively, 
and  R\  and  R2  are  the  thermal  resistances  of  the  air  and  fuel  flow, 
respectively.  Note  that  R\  and  R2  are  changing  with  respect  to 
the  air  and  fuel  flow  rates.  Therefore,  they  are  represented  as 
adjustable  resistances  in  the  equivalent  electric  circuit  model. 

The  heat  flow  Q  generated  by  chemical  reactions  can  be  cal¬ 
culated  by 


Q  =  AH- PQ 


where  AH  is  the  total  enthalpy  change  caused  by  chemical  reac¬ 
tions  and  Pe  is  the  power  consumed  by  the  electrical  loads. 

Inside  the  fuel  cell  stack,  there  are  air  and  fuel  that  have  heat 
capacities  of  Ca ir  and  Cfuei,  respectively.  Some  of  the  heat  will 
dissipate  to  the  ambient  and  the  rest  of  the  heat  goes  into  the 
SOFC  stack.  A  complete  equivalent  electric  circuit  model  that 
describes  the  whole  heat  transfer  process  is  shown  in  Fig.  2(b). 
Aggregating  the  distributed  models  together,  a  lumped  model 
is  obtained,  as  shown  in  Fig.  2(c).  R3  is  the  thermal  resistance 
of  the  heat  flow  averaged  over  the  entire  stack,  and  R4  is  the 
thermal  resistance  of  the  fuel  cell  stack. 

The  differential  equations  describing  the  lumped  model  are 
written  as 


(Cfuei  T  Cair) 


d2jnternal 
d  t 


7T  r m  t  r  m  t  r m  t  r  m  t  r  m  t 

air  1  Internal  x  *  fuel  1  Internal  x  1  sofc  1  Internal  x  ^ 

1  ~  1  ~  I 


R\ 


Csofc 


dr, 


sofc 


d  t 


r2 

^Internal  ^sofc  ^0  Csofc 


^3 


^3 


R4 


(3) 


in  which,  we  have 


Cfuei  —  ^fuel^vf  —  FfuelPfuel^vf? 


capacitances  of  the  air,  the  fuel,  and  the  SOFC  stack,  respec¬ 
tively;  raa,  raf,  and  mce  11  the  mass  of  the  air,  the  fuel,  and  one 
SOFC  cell,  respectively;  cva,  cvf,  and  cvsofc  are  specific  heat  of 
the  air,  the  fuel,  and  the  SOFC  stack,  respectively;  ns tack  the 
number  of  cells  in  a  stack;  pair  and  pfuei  the  density  of  the  air 
and  fuel,  respectively. 

Let  Ci  =  CfUei  +  Cair  and  C2  =  Csofc.  The  state  space  represen¬ 
tation  of  the  model  is  then  written  as 


X  =  AX  +  Bu,  y  =  CX  +  Du 


(4) 


X  = 


A= 


^Internal 

Csofc 

1 


X  = 


^Internal 

Csofc 


1  1 
+ 


C1C1  R2C1  R3C1 

1 


air 


U  = 


r 

Cfuei 

C0 

L  Q  J 


1 


c3Ci 


1 


1 


B  = 


1 


R3C2 

1 


R3C2  R4C2  J  . 


D  = 


C1C1  R2C1 
0  0 


0  0  0  0 
0  0  0  0 


0 

1 


1 


C4C2 


Cl 

0 


c  = 


1  0 
0  1 


If  the  SOFC  operating  temperature  is  between  650  and 
850  °C,  one  needs  to  heat  the  SOFC  stack  from  room  temperature 
to  that  temperature  range  before  operation.  In  the  heat-up  stage, 
air  flow  through  the  cathode  is  used  to  heat  the  SOFC  stack,  and 
fuel  will  not  be  fed  into  the  SOFC  stack  until  the  stack  temper¬ 
ature  exceeds  650  °C.  The  air  flow  rate,  u;a ir,  and  temperature, 
Tair,  are  control  variables  during  the  heat-up  process.  According 
to  CFD  simulation  results,  the  temperature  gradient  of  a  SOFC 
should  not  exceed  100  °C  at  any  time  to  avoid  thermal  shock. 
Therefore,  the  air  flow  temperature  should  follow  certain  rules 
such  as: 


•  When  rsofc  <  650  °C,  the  SOFC  is  running  at  heat-up  mode, 
the  following  parameters  for  the  model  shall  be  used: 

R2  =  106Js"'k_1,  e=0Js-\ 

Tair  -  T’sotc  <  I  00  °C  (5) 


•  When  Tsofe  >  650  °C,  the  SOFC  is  operating  at  running  mode, 
the  following  parameters  for  the  model  shall  be  used: 


Ri  = 


1 


Cyf 

-  Tsoic  <  50  °C 


Q  =  AH-  Pe, 


Cair  —  ^a^va  —  F^u-p^Cva?  Csofc  —  ^cell^stackTvsofc 

where  Vair  and  Vfuei  are  volumes  of  the  air  and  fuel  inside 
the  SOFC  stack,  respectively;  Cair,  Cfuei,  and  Csofc  the  thermal 


A  simulation  result  of  the  heat-up  process  (75-750  °C,  air 
flow  20gs_1)  is  shown  in  Fig.  3  The  temperature  differences 
are  controlled  to  be  within  90-100  °C  during  the  heat-up  mode. 
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Because  the  thermal  characteristics  of  different  SOFC  stacks 
may  be  different,  experimental  data  should  be  used  to  assist  in 
tuning  the  lumped  model  parameters.  Because  the  temperature 
difference  between  7intemai  and  rsofc  varies  no  more  than  10  °C 
during  the  heat-up  process,  and  in  steady  state  7internai  and  rsofc 
reach  the  same  value,  we  only  plotted  the  SOFC  stack  tempera¬ 
ture  rSOfC. 


(species  i)  in  the  reactions  can  be  written  as  [10] 


d  PV  d  n 

=  RT 


d  t 


dt’ 


d  pi  d rii 

V—  =  RT—, 
d  t  d  t 


Pi  =  ™(n?  -  <ut  +  n]) 


where  pi  is  the  partial  pressure  of  species  i  (Pa),  rii  is  the  molar 
flow  rate  of  species  i  (mol  s-1),  n]  the  total  rate  of  production  of 
species  i  (mol  s-1),  and  V  is  the  compartment  volume  (m3). 


2. 7.2.2.  Calculate  nf  If  the  SOFC  stack  current  is  I  and  there 
are  N  cells  in  a  stack,  n ]  can  be  calculated  as 


where  Q  is  the  electric  charge  transferred  during  the  reaction, 
F  (96485.3415  C  mol-1)  the  Faraday’s  constant,  and  09  is  the 
stoichiometric  coefficients  of  species  i.  Note  that  for  the  above 
chemical  reactions,  a  =  —  1  for  H2;  a  =  1  for  H2O;  and  a  =  —0.5 
for  O2. 


2.7.2.  Modeling  electrochemical  dynamics  inside  an  SOFC 
stack 

Using  hydrogen  as  fuel,  the  chemical  reactions  inside  an 
SOFC  stack  are: 

cathode  (air  electrode)  :  O2  +  4e-  — >►  20“ , 
anode  (fuel  electrode)  :  FU  +  CT  — >  H2O  +  2e-, 
cell  reaction  :  H2  +  ^C>2  —>  H2O 

The  electrochemical  dynamic  inside  an  SOFC  stack  can  be  rep¬ 
resented  in  terms  of  the  change  in  concentrations  of  the  chemical 
species  in  the  cell  reaction.  In  our  case,  the  species  of  interest 
are  FU,  O2,  and  H2O.  For  ideal  gases,  the  concentration  changes 
are  equivalent  to  the  partial  pressure  changes  in  a  fixed  volume 
container.  Therefore,  in  the  following  sections,  starting  from  the 
ideal  gas  equation,  we  derive  a  state  space  model  to  calculate 
the  partial  pressure  changes  of  the  chemical  species  and  then  use 
Nemst  equation  to  calculate  the  voltage  changes.  An  assumption 
of  the  derivations  is  that  the  cell  current  can  change  immediately 
as  long  as  the  fuel  and  oxygen  are  not  depleted. 

2. 7.2. 7.  Calculate  partial  pressures  of  chemical  species.  To 
meet  material  design  and  efficiency  requirements,  the  pressure 
inside  the  SOFC  stack  shall  not  exceed  2-3  atm.  Therefore,  it 
is  valid  to  assume  that  gases  inside  the  SOFC  stack  behave  as 
ideal  gases.  For  an  ideal  gas,  we  have 

PV  =  nRT  (7) 

where  P  is  the  absolute  pressure,  T  the  absolute  temperature,  V 
the  volume  of  the  gas,  n  the  number  of  moles,  and  R  is  the  gas 
constant  =  8.3 145 1  kJ  kmol-1  K-1 .  Then,  at  any  operating  tem¬ 
perature  T,  the  change  in  partial  pressure  p\  of  each  component 


2. 1.2. 3.  Calculate  n°ut.  When  gas  stored  under  pressure  in  a 
closed  vessel  is  discharged  to  the  atmosphere  through  a  hole 
or  other  opening,  the  gas  velocity  through  that  opening  may  be 
choked  (i.e.,  has  attained  a  maximum)  or  non-choked.  Choked 
velocity,  which  is  also  referred  to  as  sonic  velocity,  occurs  when 
the  ratio  of  the  absolute  source  pressure  Ps  to  the  absolute  down¬ 
stream  ambient  pressure  Pd  is  equal  to  or  greater  than 


Ps  _  /  K+  i\  */(*-*) 
Pd  =  (“2“ ) 


(10) 


where  K  is  the  specific  heat  ratio  of  the  discharged  gas.  For  many 
gases,  K  ranges  from  about  1.09  to  about  1.41,  and  thus  Ps/Pd 
ranges  from  1.7  to  about  1.9,  which  means  that  choked  velocity 
usually  occurs  when  the  absolute  source  vessel  pressure  is  at 
least  1.7-1. 9  times  as  high  as  the  absolute  ambient  atmospheric 
pressure  [11-16].  For  H2  and  H2O,  it  is  reasonable  to  assume 
that  the  orifice  is  choked  because  the  partial  pressures  of  FU  and 
H2O  are  low  in  atmosphere.  For  O2,  however,  depending  on  the 
reaction  rate  and  air  flow  rate,  we  may  or  may  not  treat  the  orifice 
as  choked. 

For  a  choked  orifice,  the  mass  flow  rate  wou{  can  be  calculated 
as 


CAP 


(K—l)/(2K—2) 


K+l 


(11) 


where  wout  is  the  mass  flow  rate  (g  s-1),  C  is  the  discharge  coef¬ 
ficient  (dimensionless,  usually  about  0.72),  A  is  the  discharge 
hole  area  (m2),  K  equals  cp/cv  of  the  gas,  cp  is  the  specific  heat 
at  constant  pressure,  cv  is  the  specific  heat  at  constant  volume, 
P  is  the  absolute  source  or  upstream  pressure  (Pa),  M  is  the  gas 
molecular  weight  (g  mol-1),  R  is  the  universal  gas  law  constant, 
and  T  is  the  gas  temperature  (K). 
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Because 


»r  = 


=  kf^pt 


M 


where  Xi  is  the  molar  fraction  of  species  /,  and 


Y  (  2  \(^-D/(2^-2) 


^ valve  _ 


CA\!  RT  (  K+l  J 


Vm 


The  gas  molecular  weight  M  can  be  calculated  by 


(12) 


M anode  —  ^H2 


2KVI 

Mu2  H - -r- MU20 


n 


in 

H2 


+  Xh2O^H20 


(13) 


For  example,  if  the  SOFC  operating  temperature  is  750  °C,  the 
orifice  area  is  97.96  mm2,  and  the  fuel  is  97%  (molar  fraction) 
H2  and  3%  H2O,  based  on  the  specific  heats  of  ideal  gases  in 
Table  1  [17],  the  /Talve  is  1.4298  x  10“7  and  1.4063  x  10"7  for 
H2  and  H2O,  respectively,  under  a  fuel  utilization  rate  of  0.7. 


by  H2,  Nernst  equation  can  be  written  as 


Q  rri  ^cathode 

T/  KI  1  P°2 

FNernst  —  hi 


4  F 


„ an ode 

p  o2 


(15) 


Oxygen  partial  pressure  at  the  anode  /)q ‘,dc  is  calculated  by 
equilibrium  constant  Kv  of  the  reaction  2H2  +  O2  2H2O: 


(16) 


Therefore, 


^Nernst  — 


RT 
4 F 


In 


P\i2o 


RT 
4. F 


In 


cathode  2 

Po2  Fh2  (anode) 


p\l20  (anode) 


RT 


In  K p 


For  an  isothermal  reaction,  we  have 
n  RT 

Ei  = - In  A'p 

T  4F  p 


(17) 


(IB) 


2. 1.2.4.  The  electrochemical  dynamic  model.  Based  on  Eq.  (8), 
the  state  space  representation  of  the  model  is  then  written  as 


X  =  AX  +  Bu, 


Y  =  CX  +  Du 


where 


X  = 


B  = 


A  = 


Pu2 

^cathode 

Po2 

Ph2o 


-k^WQRT 

Fan ode 

0 

0 


RT 

Fanode 

0 

0 


u  = 


n 


n 


in 

H2 

in 

o2 


n 


in 


h2o 

/ 


0 


RT 


Fcathode 

0 


0 

-kY0f&RT 

Fcathode 

0 


0 

0 

RT 


0 


0 

Fanode 
-2krRT  ~ 

Fanode 

—krRT 

Fcathode 

2krRT 


Fanode  Fanode  - 


(14) 


Because  the  operating  temperature  range  is  between  650  and 
850  °C  for  the  lumped  model,  we  use  value  at  750  °C,  which 
is  1023  K.  Then,  the  Nernst  equation  can  be  simplified  as 


FNernst  —  N 0 


RT 

2F 


t/2  \ 

ln^^) 
Pn2o  J 


(19) 


where  FNernst  is  the  open-circuit  Nernst  potential  (in  V),  £°  the 
standard  reversible  cell  potential,  pn2  the  partial  pressure  of  H2 
at  anode,  j^h2o  the  partial  pressure  of  H2O  at  anode,  po2  the 
partial  pressure  of  O2  at  cathode,  r  the  ohmic  resistance  (in  Q), 
F  the  Faraday’s  constant  (Coulomb  per  kilomole),  T  the  stack 
temperature  (in  K),  and  No  is  the  number  of  cells  in  stack. 


2.1.2. 6.  Calculate  E°.  Kv  can  be  calculated  by  the  Lewis  equa¬ 
tion: 


AGr,T  =  -RT\nKp  Kp  =  q-^.t/RT > 

q  RT  RT — AGr  t  — AGr  j 

E4=  - In  Kv  = - hi  =  - Li 

1  4F  p  4  F  RT  4  F 


(20) 


Actually,  the  change  of  Gibb’s  free  energy  AG  is  the  negative 
value  of  maximum  electric  work  W,  where  A G=  —  W=  —qAE. 
At  standard  conditions  (298.15  K  (25  °C),  1  atm), 


ArG°  =-nF  A  E° 


(21) 


"0 

0 

0 

0" 

"1 

0 

0" 

0 

0 

0 

0 

0 

1 

0 

D  = 

0 

0 

0 

0 

0 

0 

1 

*— 

_0 

0 

0 

0_ 

2. 1.2. 5.  Calculate  the  Nernst  voltage.  Nernst  equation  allows 
us  to  calculate  the  cell  potential  of  any  galvanic  cell  for  any  con¬ 
centration.  Assuming  ideal  gases,  for  fuel  cells  that  are  powered 


where  n  is  the  number  of  moles  of  electrons  generated  in  the  cell 
reactions.  From  the  definitions  of  Gibb’s  free  energy 

G  =  H  -TS  =>  ArG°  =  A tH°  -  TArS°  (22) 


where  H  is  the  enthalpy,  S  the  entropy,  and  T  is  temperature. 
Therefore, 


A rH°  -  TArS° 
nF 


(23) 
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Table  1 


Thermodynamic  data  for  the  chemical  species  in  the  SOFC  stack 


Chemical  species 

cp  (kJkg-1  K-1)  (T=750°C) 

c„  (kj  kg- 1  K-1)  (r=750°C) 

A rH°  (J  mol- 1 )  (r=750°C) 

5°  (J  K-1  mol-1)  (7’=750°C) 

h2 

14.645 

10.521 

0 

166.8831 

O2 

1.043 

0.783 

0 

244.3544 

h2o 

1.8723 

1.4108 

-247956 

233.7132 

Relevant  thermodynamic  data  for  reactants  and  products  of  the 
chemical  reactions  have  been  listed  in  Table  1  [17].  At  750  °C, 
for  2H2  +  02  -*  2H20,  we  have  £°  =  0.9915  V 


n  —  AGr  i023  -2  x  (-191326.3) 

E°l023  = - =  - 4  - -  =  0.9915  V 


4  F 


4  x  96485 


(24) 


Note:  The  value  of  ArGfjo23  is  from  the  thermodynamic  prop¬ 
erties  table  [18].  It  is  multiplied  by  2  because  there  are  2 mol 
H2O  in  the  reaction. 

The  operating  temperature  ranges  from  650  to  850  °C 
(923-1 123  K).  Using  1023  K  as  a  reference  temperature  [8],  the 
following  formula  can  be  used  to  calculate  E°  at  different  tem¬ 
peratures: 

n  0  ArS° 

A E%  =  AE%23  +  — —  (T  -  1023)  (25) 


2.7.2. 7.  Calculate  Nernst  voltage  using  partial  pressure.  Solv¬ 
ing  the  state  space  model  developed  in  previous  sections,  we 
obtain  the  voltage  dynamic  associated  with  the  changes  of  chem¬ 
ical  species. 

For  illustration  purposes,  the  fuel  flow  rate  was  set  at  0. 1  g  s~ 1 
with  a  composition  of  97%  H2  and  3%  H2O,  and  set  an  air  flow 
rate  at  0.5  g  s_1 .  Using  the  valve  constants  calculated  in  Section 
2. 1.2. 3  and  assuming  the  initial  partial  pressure  of  each  gas  is 
1  atm,  simulation  results  of  gas  partial  pressures  and  the  Nernst 
voltage  generated  across  the  stack  are  shown  in  Figs.  4  and  5, 
respectively.  The  time  constant  can  be  tens  of  seconds.  The 
steady  state  cell  voltage  is  around  1.046  V  and  drops  to  0.97  V 
when  the  cell  current  increases  from  zero  to  10  A.  Note  that  the 
partial  pressures  of  H2  and  H2O  inside  the  SOFC  stack  are  indeed 
at  least  two  times  higher  than  their  ambient  pressures.  Therefore, 
the  previous  assumption  on  the  choked  orifice  is  valid.  However, 


for  O2,  the  assumption  may  not  hold  for  lower  oxygen  flows, 
when  the  stack  current  is  high,  as  shown  in  the  i  =  10  A  case. 

2. 1.2.8.  Calculate  the  cell  voltage.  The  total  cell  voltage  can 
be  represented  as 


V  —  UNemSt  T  Vb_V  T  Vr  +  Vanode  +  Vcathode  (26) 

where  Vb_v  is  the  Butler- Volmer  over-potential  resulting  from 
charge  transfer  at  low  current,  Vr  the  voltage  drop  caused  by 
the  material  resistances,  and  Vanode  and  VCathode  are  the  over¬ 
potentials  caused  by  diffusion  in  the  cell  anode  and  cathode. 
In  the  dynamic  model,  we  assume  that  the  Vb_v  and  Vr  do  not 
change  with  respect  to  the  change  of  chemical  species.  Vanode 
and  Vcathode  are  very  small  when  the  current  is  low;  in  other 
words,  when  the  partial  pressure  of  each  gas  is  low.  When  the 
current  increases,  the  Vanode  and  Vcathode  start  to  dominate  the 
cell  output  voltage.  During  the  modeling,  the  partial  pressures  of 
each  chemical  species  are  used  to  calculate  diffusion  coefficients 
to  account  for  the  electrochemical  dynamics  in  the  cell  reaction. 

In  the  next  section,  we  will  briefly  address  the  design  con¬ 
siderations  of  an  SOFC  stack  and  then  present  our  simulation 
results  for  a  running  SOFC  stack,  obtained  by  combining  the 
thermal  dynamic  model  discussed  in  Section  2.1.1  and  the  elec¬ 
trochemical  dynamic  model  discussed  in  Section  2.1.2. 

2.2.  The  combustor  module 

The  combustor  is  mainly  used  in  a  standalone  SOFC  unit 
to  supply  extremely  high  temperature  gas  to  the  downstream 
heat  exchanger.  The  high  temperature  gas  is  used  to  heat  fuel 
and  air  inside  the  heat  exchanger  from  their  storage  temperature 
to  an  appropriate  temperature  ready  to  be  sent  to  the  SOFC 
stack.  Besides  its  main  function,  other  benefits  of  the  combustor 


0  20  40  60  80  100  120 

Time  (s) 


Fig.  4.  The  partial  pressures  of  the  H2,  O2,  and  H2O. 


Fig.  5.  The  Nernst  voltage  in  a  cell. 
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include  combusting  the  exhaust  from  the  SOFC,  recycling  extra 
heat  and  preventing  the  exhaust  of  CO. 

Three  parameters  of  the  outlet  flow  of  the  combustor  are  con¬ 
sidered:  the  outlet  flow  rate  (gs_1),  the  molar  fraction  of  the 
outlet  flow  and  the  temperature  of  the  outlet  flow.  The  combus¬ 
tor  model  is  explained  in  two  parts  in  this  section.  The  first  part  is 
the  derivation  of  the  outlet  flow  rate  and  its  molar  fraction.  The 
second  part  is  the  calculation  of  the  temperature  of  the  outlet 
flow. 


Table  2 


Specific  heat  and  enthalpy  of  the  chemical  species  in  the  SOFC  stack 


Chemical  species 

Enthalpy  of  formation 

A  H°  (J  mol-1) 

Cp  (J  mol-1  K-1) 

h2 

0 

28.8 

o2 

0 

29.4 

h20(1) 

-285.830 

75.3 

H20(g) 

-241.820 

33.6 

n2 

0 

29.1 

2.2.7.  The  outlet  flow  rate  and  its  molar  fraction 

Assume  that  we  know  the  inlet  flow  rate  wm  (gs_1),  inlet 
flow  molar  fraction  (yg ,  yjlj  Q,  ylQ2,  yj^2),  and  inlet  temperature 

Tn  (°C).  The  molar  flow  rate  of  the  inlet  flow  can  be  calculated 
by 


nm  = 


UJ 


in 


y'&M h2  +  yli0MH2o  +  y'£M o2  +  ySMN. 


(27) 


12  ^  fl2U  '  •/U2 

The  molar  flow  rate  of  each  gas  component  is  given  by 
nf1  =  nmym 


(28) 


The  reaction  that  happened  in  the  combustor  can  be  expressed  as 
2H2  +  O2  — >  2H2O.  Because  the  SOFC  stack  and  the  combustor 
are  always  supplied  with  surplus  air,  it  is  reasonable  to  assume 
that  100%  of  the  inlet  hydrogen  is  combusted.  Then  we  have  the 
molar  changes  of  each  of  the  species  after  the  reaction: 


combusted  _  in 

wH2  -  "h2 


^consumed  =  Q  5„m 


n 


generated  _  in 


h2o  —  nn2 
The  outlet  molar  flow  rate  of  each  individual  gas  is 


(29) 


„°ut  _ 
wH2  - 

nin  — 

"h2 

M  combusted 

«h2 

0, 

Mout 

«o2  = 

Min 

«o2  - 

„  consumed 

«o2 

T/,out 

^in 

.  generated 

^in  i  Min 

wh2o  - 

-  wh2o 

+  «h2o 

—  nH2Q  +  nH2’ 

Mout 

™in 

«n2  = 

"n2 

(30) 


heat  energy  of  the  inlet  gas  and  the  energy  released  from  the 
combustion  of  hydrogen.  The  energy  balance  of  gases  inside  the 
combustor  is  then  represented  by 

gout  =  gin  +  gcom  (34) 


where  gcom  is  the  energy  of  combustion,  Qou{  the  energy  to  heat 
the  outlet  gas  from  298  K  to  the  outlet  temperature,  and  Qm  is 
the  energy  to  cool  the  inlet  gas  from  the  inlet  temperature  to 
298  K. 

The  energy  needed  to  cool  1  mol  of  the  inlet  flow  q  to  298  K 
can  be  calculated  by 

q'u2  =  cf(Tm  -  298),  =  cf(Tm  -  298), 

7h2o  =  Cp2°®(373  -  298)  +  (A//f^H20  g)  —  A//^H20  ]}) 

+  c^20(g)(rin  -  373), 

q™2=cf(Tm-  298)  (35) 

Note  that  the  term  (A/7f^H90  ^  —  A/7f^H20  ^),  which  equals 

44  010  J  mol-1,  appears  in  the  equation  as  the  energy  needed 
to  gasify  the  water,  or  the  molar  enthalpy  of  vaporization  of 
water.  Refer  to  Table  2  for  the  specific  heat  and  enthalpy  of  the 
chemical  species  in  an  SOFC  stack. 

The  energy  carried  by  the  inlet  flow  is  given  by 


(36) 


The  total  outlet  molar  flow  rate  is  the  summation  of  the  molar 
flow  rate  of  individual  gases: 


„  out  „  out  I  „  out  I  „  out  ,  „  out 

n  =  «H2  +  no2  +  wh2o  +  «n2 


02  T  "H20  t  “N2  (31) 

The  molar  fraction  of  the  outlet  flow  can  be  calculated  by 


yfut  = 


nr 


n 


out 


(32) 


Finally  we  get  the  mass  flow  rate  (g  s  l)  of  the  outlet  flow: 


w 


out 


=  <Mh2  +  n(Q2Mo2  +  ^H2L0MH20  +  /i^Mno  (33) 


out 


.out 


out 


2.2.2.  Temperature  calculation  of  the  outlet  flow 

Temperature  calculation  of  the  outlet  flow  is  based  on  energy 
balance.  The  energy  released  from  combustion  is  equal  to  the 
energy  absorbed  by  the  outlet  gases.  Using  25  °C  (298  K)  as  a 
reference  for  potential  heat  energy,  the  potential  heat  energy  of 
the  outlet  gas  should  be  equal  to  the  summation  of  the  potential 


The  energy  Q  released  by  the  combustion 


gcom 


gen 

h2o 


AHm2oA) 


—  n 


com 

H2 


A^H2) 


(37) 


The  energy  carried  by  1  mol  of  the  outlet  flow  can  be  calculated 
by 

qf  =  cf(Tout  -  298),  q'f  =  c°2(rout  -  298), 

=  c^2°(1)(373  -  298)  +  44010  +  c^CKg^out  _  373 ^ 

=  f,2(Tout  -  298)  (38) 

The  energy  carried  by  the  outlet  flow  is  given  by 

<2out  =  «  +  «  +  <o<o  +  «n2‘C  (39) 


T°ut  can  then  be  calculated  by  solving  the  above  equations. 
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2.3.  The  heat  exchanger  module 

A  heat  exchanger  is  used  in  the  SOFC  unit  to  heat  the  inflow 
gases  to  an  appropriate  temperature  using  the  heat  generated 
by  combustion  of  the  inflow  gases  and  the  heat  carried  by  the 
exhaust  gases.  Heat  is  transferred  from  the  hot  inlet  flow  to 
the  cold  inlet  flow.  There  are  two  major  concerns  here  in  the 
modeling  of  a  heat  exchanger:  the  heat  exchanging  rate  and  the 
maximum  heat  that  can  be  transferred.  In  the  following  discus¬ 
sion,  annotation  4 1  ’  stands  for  the  hot  flow,  whereas  annotation 
‘2’  represents  the  cold  flow  to  be  heated.  The  heat  exchange 
process  can  be  represented  by 

=  Gin  -  QT  ~  Gtransfer, 

PiCpi  V2  2  r20Ut  =  Qf  -  QT  +  G,ransfer  (40) 

So  we  have 

PlcplVijTr  =  wiCpi(r;n  -  tT)  -  Gtransfer, 

P2Cp2V2jTr  =  w2Cp2(T f  -  7T)  +  Gtransfer  (41) 

where  Wi  is  the  mass  flow  rate  (kgs-1),  cPi  the  specific  heat 
(Jkg-1  K-1),  and  Qi  is  the  heat. 

The  heat  transfer  effectiveness  is  defined  as  the  actual  heat 
transfer  rate  divided  by  the  maximum  possible  heat  transfer  rate 

[19]: 

q  transfer 

e  =  —. -  (42) 

Q  max 

Let  us  assume  the  heat  transferred  between  the  two  gases  is  a 
function  of  the  inlet  temperatures  in  the  following  form: 

^transfer  =  ^in  _  rm}  (43) 

Then  the  heat  exchange  rate  rises  to  its  maximum  value  when 
kq  takes  its  maximum  value 

Qmax  =  kfax\T[n  -  if)  (44) 

The  maximum  heat  exchange  rate  happens  when  the  outlet  tem¬ 
peratures  of  the  two  gases  are  equal.  Solving  the  differential 
equations,  the  outlet  temperatures  can  be  expressed  as  functions 
of  time 


Fig.  6.  An  example  of  the  calculation  results  of  the  heat  exchanger  model. 


When  7jout  equals  T™1,  the  maximum  heat  exchange  rate  hap¬ 
pens  and  kq  takes  its  maximum  value: 


U)\C p\  W2Cp2 
WiCpi  +  W2Cp2 


So, 

_  wlcplw2Cp2  prpin 
kfmax  —  ,  f 1 1 

W\Cpi  +  W2Cp2 


Let  Cmjn  — 


We  get 


W\Cp\  W2Cp2 
WiCpi+W2Cp2  ’ 


P\cp\V\jT"m  =  w\cp\(Tm 
P2Cp2V2jT^  =  W2Cp2(T/ 


77')  -  sCmin(T/ 


r20Ut)  +  sCmin(T[ 


in 


(47) 


(48) 


Fig.  6  shows  an  example  of  the  calculation  results  of  the  heat 
exchanger  model  when  the  heat  transfer  coefficient  is  equal 
to  0.7  and  1.0,  respectively.  The  temperature  of  the  hot  inlet 
flow  is  1 300  K  in  the  example,  whereas  the  temperature  of  the 
cold  inlet  flow  is  315  K.  The  calculation  results  show  that  when 
heat  transfer  coefficient  is  equal  to  1 .0,  the  outlet  temperatures 


W\Cp\  kq  ^ in  ^ 

W\Cp\  1 

U)2Cp2  ~  kq  ^in  | 
l2  "h 

U)2CP2 


kq 

W\Cp\ 


kq 

W2Cp2 


W\Cp\  kq  t in 

W\Cp\  1 

U)2Cp2  ~  kq  ^in  | 
l2  "b 

0J2CP2 


kq 

W\Cp\ 


kq 

W2Cp2 


Q(-wicpi/picpiVi)t 

Q(-W2Cp2/p2Cp2V2)t 


The  steady  state  solution  of  the  differential  equation  is 


(45) 


T^OUt 

71 


Ql\Cp\  kq  t in  ^ 

WiCpi  1 


kq 

W\Cpi 


T-'in 
l2  ’ 


0)2C  p2  kq 


kq 


(46) 


of  both  flows  are  equal.Based  on  the  above  equations,  the 
state  space  representation  of  the  partial  pressure  model  is  then 
written  as 


W2CP2 


W2CP2 


X  =  AX  +  Bu, 


Y  =  CX  + Du 


(50) 
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2.4.  The  simulation  results  of  an  SOFC  APU  at  running 
mode 
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In  the  controller  design  for  an  SOFC  stack,  there  should  be 
three  control  objectives:  temperature,  voltage  and  fuel  utiliza¬ 
tion.  The  following  constraints  have  to  be  met: 

(A)  The  stack  temperature  will  not  exceed  its  upper  and  lower 
temperature  limits. 

(B)  The  voltage  should  be  managed  to  be  above  its  lower  limit 
while  A  is  satisfied. 

(C)  The  fuel  utilization  should  be  above  its  lower  limit  provided 
that  (A)  and  (B)  are  satisfied. 

Because  it  is  not  the  focus  of  this  paper  to  discuss  the  control 
scheme  of  the  SOFC  APU,  a  simplified  controller  model  is  cho¬ 
sen.  It  consists  of  two  parts:  the  control  of  the  fuel  utilization  and 
the  control  of  the  heat  exchanger  flow  mix.  With  the  controller, 
one  can  control  the  fuel  and  air  flow  rates  and  temperatures. 

Using  MATLAB,  we  implemented  the  above  thermal 
dynamic  and  electrochemical  dynamic  models  to  simulate  an 
SOFC  stack  containing  60  cells  at  operation  mode.  The  load 
is  purely  resistive.  At  5  s,  the  load  changes  from  60  to  3  Q.  As 
shown  in  Fig.  7,  the  voltage  then  drops  from  60  to  41.3  V.  If 
we  set  the  fuel  flow  to  ramp  up  at  0.04  g  s-1  from  15  to  20  s,  as 
shown  in  Fig.  8,  we  have  the  voltage  and  current  increase  accord¬ 
ingly.  However,  the  fuel  utilization  will  drop  significantly  from 
0.55  down  to  0.06. 

From  the  simulation  results  shown  in  Figs.  7  and  8,  we  have 
the  following  observations: 

•  Different  dynamics  have  different  time  constants.  As  shown 
in  Fig.  7,  the  stack  temperature  will  not  change  significantly 
during  tens  of  seconds;  therefore,  thermal  dynamics  will  not  in 
general  affect  the  electric  dynamic  responses,  which  happen 
in  a  few  seconds  at  the  most.  Electrochemical  dynamics  are 
within  seconds  and  may  interact  with  the  electric  dynamics 
caused  by  outside  circuits.  The  fluid  dynamics  may  happen 
in  seconds  because  of  the  fan,  valves,  and  the  inlet  and  outlet 
parameters. 

•  When  fuel  flow  rates  increase,  one  can  boost  the  stack  volt¬ 
age  and  current  in  seconds.  However,  the  fuel  utilization  will 


Fig.  7.  The  stack  temperatures,  voltage,  and  current. 
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Fig.  8.  The  fuel  flow  rate,  the  fuel  utilization,  and  partial  pressures  of  an  SOFC 
stack. 


drop  sharply,  and  the  pressure  inside  the  SOFC  stack  will 
increase  as  well,  causing  low  efficiency  and  high  mechanical 
stresses  inside  the  SOFC,  both  of  which  are  not  desired.  Fur¬ 
thermore,  the  rise  in  voltage  is  not  significant  when  fuel  flow 
rate  increases  further.  Therefore,  the  voltage  regulation  of  the 
SOFC  is  very  limited. 

3.  Conclusions 

This  paper  has  addressed  the  dynamic  modeling  issues  in 
modeling  a  standalone  SOFC  auxiliary  power  unit.  While  pro¬ 
viding  greater  accuracy,  detailed  SOFC  models  such  as  CFD 
take  considerable  computer  resources  to  run,  and  it  is  inconve¬ 
nient  to  couple  detailed  models  with  an  electric  circuit  model 
to  simulate  the  electric  circuit  dynamic  response.  An  objective 
of  this  research  has  been  to  set  up  lumped  thermal  dynamic  and 
electrochemical  dynamic  models  to  improve  the  computational 
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time  and  provide  reasonable  accuracy  of  the  calculated  stack 
voltage,  temperature,  and  fuel  utilization.  Detailed  derivation 
of  model  parameters  has  been  presented.  The  constraints  and 
parameters  of  the  lumped  model  have  been  obtained  through 
detailed  models  and  experimental  data.  This  SOFC  dynamic 
model  is  expected  to  be  used  for  SOFC  controller  design  as 
well  as  simulating  the  dynamic  response  of  SOFC  in  power 
systems. 
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